***************************************************************
* .do-file to estimate income process of individuals from PSID  
***************************************************************

clear all
set maxvar 7000

use "$psidwkddata/famind_combined_adj_sing.dta", clear

///////////////////////////
/////// HOUSEKEEPING //////
///////////////////////////

// generate personal identifier
gen uniqueid_ref = (famid68*1000) + persid if ref    == 1
gen uniqueid_sp  = (famid68*1000) + persid if spouse == 1

gen uniqueid     = uniqueid_ref if ref == 1
replace uniqueid = uniqueid_sp  if spouse == 1

//////////////////////////////////////
/////// RECODE INCOME VARIABLES //////
//////////////////////////////////////

sum inc*
replace inc_rent_ref      = 0 if inc_rent_ref      == -99999 | inc_rent_ref      == 999998 | inc_rent_ref == 999999
replace inc_div_ref       = 0 if inc_div_ref       == 999998 | inc_div_ref       == 999999
replace inc_interest_ref  = 0 if inc_interest_ref  == 999998 | inc_interest_ref  == 999999
replace inc_trustfund_ref = 0 if inc_trustfund_ref == 999998 | inc_trustfund_ref == 999999

// re-label income variables from missing to 0 
*total 
replace totfaminc  	= 0 if totfaminc == .
replace inc_ss_all  	= 0 if inc_ss_all == .

*other family members
replace inc_taxable_oth = 0 if inc_taxable_oth == .
replace inc_asset_oth   = 0 if inc_asset_oth == .
replace inc_labor_oth   = 0 if inc_labor_oth == .
replace inc_transf_oth  = 0 if inc_transf_oth == .
replace inc_ss_oth      = 0 if inc_ss_oth == .

*reference person + spouse
replace inc_taxable_refsp = 0 if inc_taxable_refsp == .
replace inc_transf_refsp  = 0 if inc_transf_refsp == .

*reference person
replace inc_ss_ref		  = 0 if inc_ss_ref   == .
replace inc_rent_ref      = 0 if inc_rent_ref == . 
replace inc_div_ref       = 0 if inc_div_ref  == . 
replace inc_interest_ref  = 0 if inc_interest_ref == . 
replace inc_trustfund_ref = 0 if inc_trustfund_ref == . 
replace inc_labor_ref     = 0 if inc_labor_ref == .
replace inc_laborbus_ref  = 0 if inc_laborbus_ref == . 

*spouse
replace inc_ss_sp         = 0 if inc_ss_sp == .
replace inc_rent_sp       = 0 if inc_rent_sp == . 
replace inc_div_sp        = 0 if inc_div_sp == . 
replace inc_interest_sp   = 0 if inc_interest_sp == . 
replace inc_trustfund_sp  = 0 if inc_trustfund_sp == .
replace inc_labor_sp      = 0 if inc_labor_sp == .
replace inc_laborbus_sp   = 0 if inc_laborbus_sp == . 

////////////////////////////////////////////////////////////////////////////

// HH income (excluding asset income), sep. for reference person and spouse

////////////////////////////////////////////////////////////////////////////

gen inc_noasset_ref = 0
gen inc_noasset_sp  = 0

foreach num of numlist 1989/1993 {
// couples
replace inc_noasset_ref = inc_labor_ref + inc_laborbus_ref + inc_ss_ref + ((inc_ss_oth + inc_transf_oth + inc_transf_refsp + inc_labor_oth)/2)  ///
			if wave == `num' & single == . 
			
			
replace inc_noasset_sp =  inc_labor_sp + inc_laborbus_sp +  inc_ss_sp + ((inc_ss_oth + inc_transf_oth + inc_transf_refsp + inc_labor_oth)/2) ///
			if wave == `num' & single == . 

// singles
replace inc_noasset_ref = inc_labor_ref + inc_laborbus_ref + inc_ss_ref + inc_ss_oth + inc_transf_oth + inc_transf_refsp + inc_labor_oth   ///
			if wave == `num' & single != . 
}


foreach num of numlist 2005(2)2017 {
// couples
replace inc_noasset_ref = inc_labor_ref + inc_laborbus_ref + inc_ss_ref + ((inc_ss_oth + inc_transf_oth + inc_transf_refsp + inc_labor_oth)/2)  ///
			if wave == `num' & single == . 
			
			
replace inc_noasset_sp =  inc_labor_sp + inc_laborbus_sp + inc_ss_sp + ((inc_ss_oth + inc_transf_oth + inc_transf_refsp + inc_labor_oth)/2) ///
			if wave == `num' & single == . 

// singles
replace inc_noasset_ref = inc_labor_ref + inc_laborbus_ref + inc_ss_ref + inc_ss_oth + inc_transf_oth + inc_transf_refsp + inc_labor_oth   ///
			if wave == `num' & single != . 
}



foreach num of numlist 1994/1996 {
// couples
replace inc_noasset_ref =  inc_labor_ref + inc_laborbus_ref + ((inc_ss_all + inc_transf_oth + inc_transf_refsp + inc_labor_oth)/2)  ///
							if wave == `num' & single == .
							
replace inc_noasset_sp  =  inc_labor_sp + inc_laborbus_sp   + ((inc_ss_all + inc_transf_oth + inc_transf_refsp + inc_labor_oth)/2)  ///
							if wave == `num' & single == .

// singles
replace inc_noasset_ref =  inc_labor_ref + inc_laborbus_ref + inc_ss_all + inc_transf_oth + inc_transf_refsp + inc_labor_oth  ///
							if wave == `num' & single != .
}


foreach num of numlist 1997(2)2003 {
// couples
replace inc_noasset_ref =  inc_labor_ref + inc_laborbus_ref + ((inc_ss_all + inc_transf_oth + inc_transf_refsp + inc_labor_oth)/2)  ///
							if wave == `num' & single == .
							
replace inc_noasset_sp  =  inc_labor_sp + inc_laborbus_sp + ((inc_ss_all + inc_transf_oth + inc_transf_refsp + inc_labor_oth)/2)  ///
							if wave == `num' & single == .

// singles
replace inc_noasset_ref =  inc_labor_ref + inc_laborbus_ref + inc_ss_all + inc_transf_oth + inc_transf_refsp + inc_labor_oth  ///
							if wave == `num' & single != .
}


//////////////////////////////////////
/////// FURTHER DATA CLEANING ////////
//////////////////////////////////////

// recode sexind as dummy 
replace sexind = 0 if sexind == 1
replace sexind = 1 if sexind == 2

// winsorize (as for asset data) 
winsor2 inc_noasset_ref, replace cuts(1 99)
winsor2 inc_noasset_sp, replace cuts(1 99)

// drop all observations that are neither reference person nor spouse
keep if uniqueid != .

// create panel component
xtset uniqueid wave

// indvidual income variable 
gen incind = 0
replace incind = inc_noasset_ref if ref == 1
replace incind = inc_noasset_sp  if spouse == 1

// create additional variables for regression (all from individual perspective)
gen genderage = sexind*ageind
gen ageindsq  = ageind^2
gen ageindcb  = ageind^3

// restrict to sample to individual age 30 - 65
keep if ageind <= 65 & ageind >= 30 

/////// WHEN ESTIMATING DET. INCOME //////////
// drop all single individuals who report zero income and couples who report zero total income
replace incind = . if inc_noasset_ref == 0 & inc_noasset_sp == 0 & single == . 
replace incind = . if incind == 0 & single !=.

// transformation of income variables
ihstrans incind   // work with ihs bc of small income values

// constant weights within individuals
bys uniqueid: egen panelweight=mean(famweight)

//////////////////////////////////////////////////////////////////////////////////////////
//////////////////////////// LIFE CYCLE PROFILES OF DET INCOME ///////////////////////////
//////////////////////////////////////////////////////////////////////////////////////////

// COUPLES 
xtreg ihs_incind ageind ageindsq genderage if single == . & educcat != . [aweight = panelweight] ,fe vce(robust)	

predict pr 	  if ihs_incind != . & single == . & educcat != .
predict resid if ihs_incind != . & single == . & educcat != ., ue    
bysort ageind sexind: egen mean = mean(pr) if single == .
drop pr

outreg2 using "$resultpath/detearn_coup.tex", ///
ctitle("C -- first stage")  replace

reg resid i.educcat##sexind if ihs_incind != . & single == . [aweight = panelweight], vce(robust) 
predict pr0 if ihs_incind != . & single == .	

outreg2 using "$resultpath/detearn_coup.tex", ///
ctitle("C -- second stage") append

bysort ageind sexind: egen mean0 = mean(pr0) if single == .	

gen mean1     = mean + mean0      if ihs_incind != . & single == .	
replace mean1 = exp(mean1-log(2)) if ihs_incind != . & single == .	
replace mean1 = mean1/1000        if ihs_incind != . & single == .		

sort ageind
tw  (line mean1 ageind   if sexind == 0 & single == ., lc(black) lw(medthick)) ///
    (line mean1 ageind   if sexind == 1 & single == ., lp(dash) lc(black) lw(medthick)),  ///
xtitle("Age", size(large)) ytitle("(000s) annual income in 2007 \$", size(large) height(10))  ///
legend(order(1 "Men" 2 "Women") size(vlarge) position(6) col(4) region(lcolor(white)) stack) ///		
graphregion(color(white)) xlabel(30(10)65,labsize(vlarge)) ylabel(,labsize(vlarge))
graph export "$resultpath/detincome_couples.eps", replace
graph export "$resultpath/detincome_couples.pdf", replace

drop mean resid  pr0 mean* 


/// SINGLES (single men and single women)
xtreg ihs_incind ageind ageindsq genderage if single != . & educcat != . [aweight = panelweight],fe vce(robust)	
		
predict pr if ihs_incind     != . & single != .	& educcat != .
predict resid if ihs_incind  != . & single != . & educcat != ., ue    
bysort ageind sexind: egen mean = mean(pr) if single != .	
drop pr

outreg2 using "$resultpath/detearn_sing.tex", ///
ctitle("C -- first stage")  replace

reg resid i.educcat##sexind if single != .	[aweight = panelweight], vce(robust) 
predict pr0 if ihs_incind   != . & single != .	

outreg2 using "$resultpath/detearn_sing.tex", ///
ctitle("C -- second stage") append

bysort ageind sexind: egen mean0 = mean(pr0) if single != .
	
gen mean1   = mean + mean0   if single != .	
replace mean1 = exp(mean1-log(2))  if single != .	
replace mean1 = mean1/1000     if single != .

sort ageind
tw  (line mean1 ageind   if sexind == 0 & single != ., lc(black) lw(medthick)) ///
    (line mean1 ageind   if sexind == 1 & single != ., lp(dash) lc(black) lw(medthick)),  ///
xtitle("Age", size(large)) ytitle("(000s) annual income in 2007 \$", size(large) height(10))  ///
legend(order(1 "Men" 2 "Women") size(vlarge) ///
		position(6) col(4) region(lcolor(white)) stack) ///
graphregion(color(white)) xlabel(30(10)65,labsize(vlarge)) ylabel(,labsize(vlarge))
graph export "$resultpath/detincome_singles.eps", replace
graph export "$resultpath/detincome_singles.pdf", replace

drop mean resid  pr0 mean* 

//////////////////////////////////////////////////////////////////////////////////////////
//////////////////////////// TRANSITORY INCOME COMPONENT /////////////////////////////////
//////////////////////////////////////////////////////////////////////////////////////////

drop ihs_incind

////////// WHEN ESTIMATING STOCH. INCOME /////////////
// drop all individuals who report zero income: done to avoid very large variances for married women (who drop out and re-enter lf)
// note: does not affect singles
replace incind = . if incind == 0

// transformation of income variables
ihstrans incind   // work with ihs bc of small income values

// COUPLES 
xtreg ihs_incind ageind ageindsq genderage if single == . & educcat != . [aweight = panelweight] ,fe vce(robust)	

predict pr if ihs_incind != . & single == .	& educcat != .
predict resid if ihs_incind  != . & single == . & educcat != ., ue    
bysort ageind sexind: egen mean = mean(pr) if single == .
drop pr

reg resid i.educcat##sexind if ihs_incind != . & single == . [aweight = panelweight], vce(robust) 
predict pr0 if ihs_incind != . & single == .	

gen resid11  = resid-pr0  if sexind == 0  & ihs_incind != . & single == .	
gen resid12  = resid-pr0  if sexind == 1  & ihs_incind != . & single == .	

drop mean resid  pr0 mean* 


/// SINGLES 
xtreg ihs_incind ageind ageindsq genderage if single != . & educcat != . [aweight = panelweight],fe vce(robust)	
		
predict pr if ihs_incind     != . & single != .	& educcat != .
predict resid if ihs_incind  != . & single != . & educcat != ., ue    
bysort ageind sexind: egen mean = mean(pr) if single != .	
drop pr

reg resid i.educcat##sexind if single != .	[aweight = panelweight], vce(robust) 
predict pr0 if ihs_incind   != . & single != .	

gen resid13  = resid -pr0  if sexind == 0  & single != .	
gen resid14 = resid -pr0  if sexind == 1   & single != .
	
drop mean resid  pr0 mean* 


*********** EXPORT RESIDUALS TO ESTIMATE INCOME PROCESS IN MATLAB ********

// married men
preserve
keep if sexind == 0 & single == .
keep wave uniqueid resid11 

reshape wide resid11, i(uniqueid) j(wave)
outsheet  using "$resultpath/COVmale_coup_labinc.csv",comma nolabel  nonames replace
restore

// married women
preserve
keep if sexind == 1 & single == . 
keep wave uniqueid resid12

reshape wide resid12, i(uniqueid) j(wave)
outsheet  using "$resultpath/COVfemale_coup_labinc.csv",comma nolabel  nonames replace
restore


// single men 
preserve
keep if sexind == 0 & single != .
keep wave uniqueid resid13 

reshape wide resid13, i(uniqueid) j(wave)
outsheet  using "$resultpath/COVmale_sing_labinc.csv",comma nolabel  nonames replace
restore

// single women 
preserve
keep if sexind == 1 & single != . 
keep wave uniqueid resid14

reshape wide resid14, i(uniqueid) j(wave)
outsheet  using "$resultpath/COVfemale_sing_labinc.csv",comma nolabel  nonames replace
restore

preserve
keep uniqueid wave resid*

save "$resultpath/incomeresid.dta", replace
restore


